EFECTO DE LA CALIDAD DE AIRE EN DISTINTAS ENFERMEDADES.

1.Introducción

En este seminario estamos en búsqueda de la relación entre los niveles de diferentes parámetros de los componentes del aire y las muertes que ocurren en España por comunidades autónomas diferenciadas por su causa.

El resultado que pretendemos obtener es que cuanto mayor sean los niveles de los diferentes parámetros del aire mayor deberá ser el numero de muerte, esto en principio indicara que hay una cierta relación en que la muerte puede estar causada por una enfermedad que ha sido expresada por un nivel alto de PM 2.5 o arsénico, por ejemplo.

De las 4 preguntas que hemos realizado, previos estudios nos indican que hay una cierta correlación entre los niveles de arsénico y el cáncer, esto ocurre porque el arsénico su principal efecto nocivo será el cáncer de piel, esto veremos si se ve reflejado en nuestro estudio. También ocurre con el estudio de como se ve afectada la neumonía por el PM 2,5 ya que en estudios testados nos indican que una mayor exposición a contaminantes ambientales como el PM 2,5 y PM10 es un factor de riesgo para pacientes por neumonía.

Los datos los hemos obtenido a partir de las páginas web puestas a continuación:

2.Objetivos

El objetivo principal es encontrar cualquier efecto de la calidad de aire entre las distintas enfermedades:

  • Encontrar el efecto de la cantidad de PM 2,5 sobre las muertes por Neumonía.

  • Encontrar el efecto de la cantidad de PM 10 sobre las muertes por 4 enfermedades, Asma, Gripe, Neumonia y Enfermedades en las Vias Inferiores.

  • Encontrar el efecto de la cantidad de Arsénico sobre las muertes por diabetes.

  • Encontrar algún tipo de correlación entre el Benzopireno y Arsénico sobre las muertes producidas por Cáncer. Ver cual de las dos variables tiene más influencia sobre las muertes producidas por Cáncer.

3.Importar librerias

Comenzamos importando todas las librerías con las que vamos a trabajar en el documento.

library(readr)
library(dplyr)
library(tidyverse)
library(tidyselect)
library(ggplot2)
library(DT)

4.Datos de Calidad de aire.

En este caso vamos a importar los datos de Calidad del aire de forma separada, se encuentran en varias tablas dependiendo del año y las magnitudes de calidad.

Año 2020

Comenzamos con el año 2020:

#2020 ---------------------------------------------------------------------------------------------------
As_DD_2020 <- read_delim("input/data/calidad del aire/2020/As_DD_2020.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

BaP_DD_2020 <- read_delim("input/data/calidad del aire/2020/BaP_DD_2020.csv", 
                          delim = ";", escape_double = FALSE, trim_ws = TRUE)

Cd_DD_2020 <- read_delim("input/data/calidad del aire/2020/Cd_DD_2020.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Ni_DD_2020 <- read_delim("input/data/calidad del aire/2020/Ni_DD_2020.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Pb_DD_2020 <- read_delim("input/data/calidad del aire/2020/Pb_DD_2020.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM10_DD_2020 <- read_delim("input/data/calidad del aire/2020/PM10_DD_2020.csv", 
                           delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM25_DD_2020 <- read_delim("input/data/calidad del aire/2020/PM25_DD_2020.csv", 
                           delim = ";", escape_double = FALSE, trim_ws = TRUE)

DT::datatable(As_DD_2020)

Ahora tenemos que modificar las distintas tablas importadas anteriormente para conseguir tratar los datos de calidad de aire en una única tabla.

##Arsenico 2020
tablaAS_2020 <- select(.data = As_DD_2020, PROVINCIA, ANNO, D01:D31 )

tablaAS_2020 <- mutate(tablaAS_2020, AÑO = 2020)

tabAS_2020 <- tablaAS_2020 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabAS_2020$ARSENICO <- apply(tabAS_2020[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabAS_2020 <- select(.data = tabAS_2020, PROVINCIA, AÑO, ARSENICO)

##Benzopireno 2020
tablaBaP_2020 <- select(.data = BaP_DD_2020, PROVINCIA, ANNO, D01:D31)

tablaBaP_2020 <- mutate(tablaBaP_2020, AÑO = 2020)

tabBaP_2020 <- tablaBaP_2020 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabBaP_2020$BENZOPIRENO <- apply(tabBaP_2020[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabBaP_2020 <- select(.data = tabBaP_2020, PROVINCIA, AÑO, BENZOPIRENO)

##Cadmio 2020

tablaCd_2020 <- select(.data = Cd_DD_2020, PROVINCIA, ANNO, D01:D31)

tablaCd_2020 <- mutate(tablaCd_2020, AÑO = 2020)

tabCd_2020 <- tablaCd_2020 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabCd_2020$CADMIO <- apply(tabCd_2020[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabCd_2020 <- select(.data = tabCd_2020, PROVINCIA, AÑO, CADMIO)

##Niquel 2020

tablaNi_2020 <- select(.data = Ni_DD_2020, PROVINCIA, ANNO, D01:D31)

tablaNi_2020 <- mutate(tablaNi_2020, AÑO = 2020)

tabNi_2020 <- tablaNi_2020 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabNi_2020$NIQUEL <- apply(tabNi_2020[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabNi_2020 <- select(.data = tabNi_2020, PROVINCIA, AÑO, NIQUEL)

##Plomo 2020

tablaPb_2020 <- select(.data = Pb_DD_2020, PROVINCIA, ANNO, D01:D31)

tablaPb_2020 <- mutate(tablaPb_2020, AÑO = 2020)

tabPb_2020 <- tablaPb_2020 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPb_2020$PLOMO <- apply(tabPb_2020[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPb_2020 <- select(.data = tabPb_2020, PROVINCIA, AÑO, PLOMO)

##PM10 2020

tablaPM10_2020 <- select(.data = PM10_DD_2020, PROVINCIA, ANNO, D01:D31)

tablaPM10_2020 <- mutate(tablaPM10_2020, AÑO = 2020)

tabPM10_2020 <- tablaPM10_2020 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM10_2020$PM_10 <- apply(tabPM10_2020[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM10_2020 <- select(.data = tabPM10_2020, PROVINCIA, AÑO, PM_10)

##PM25 2020

tablaPM25_2020 <- select(.data = PM25_DD_2020, PROVINCIA, ANNO, D01:D31)

tablaPM25_2020 <- mutate(tablaPM25_2020, AÑO = 2020)

tabPM25_2020 <- tablaPM25_2020 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM25_2020$PM_25 <- apply(tabPM25_2020[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM25_2020 <- select(.data = tabPM25_2020, PROVINCIA, AÑO, PM_25)

DT::datatable(tabAS_2020)

Por cada magnitud de calidad de aire obtenemos una tabla con tres colummnas: provincia, año y dato. De la siguiente manera.En caso del año 2020 hemos tenido que realizar una modificación que no hemos hecho en las demas para cambiar el año. Por cada magnitud, lo que hacemos es seleccionar de la tabla original las columnas que nos hacen falta y creando otra tabla a partir de estas. Después, creamos otra tabla haciendo un group by y un summarise para hacer la media de los dias dependiendo del año y de las provincias. Por último, hacemos la media de todos los dias de un mes por cada provincia y año. Obtenemos una colunmna con los datos de nuestra magnitud por años y por provincia.

Union de los datos del año 2020.

Para unir las 7 tablas en una única realizamos un full_join de dos en dos. De la siguiente manera:

#union de los datos
join2020.1 <- full_join (x=tabAS_2020,y=tabBaP_2020, by=c("PROVINCIA","AÑO"))

join2020.2 <- full_join (x=join2020.1,y=tabCd_2020, by=c("PROVINCIA","AÑO"))

join2020.3 <- full_join (x=join2020.2,y=tabNi_2020, by=c("PROVINCIA","AÑO"))

join2020.4 <- full_join (x=join2020.3,y=tabPb_2020, by=c("PROVINCIA","AÑO"))

join2020.5 <- full_join (x=join2020.4,y=tabPM10_2020, by=c("PROVINCIA","AÑO"))

join2020.6 <- full_join (x=join2020.5,y=tabPM25_2020, by=c("PROVINCIA","AÑO"))

DT::datatable(join2020.6)

A partir del código anterior obtenemos una tabla con 9 columnas que se corresponden a: PROVINCIA, AÑO, ARSENICO, BENZOPIRENO, CADMIO, NIQUEL, PLOMO, PM_10, PM_25.

Años 2019, 2018, 2017 y 2016

Con los otros 4 años realizamos lo mismo. Lo único que cambia es que a la columna “ANNO” de la tabla original le cambiamos el nombre a “AÑO” para evitar errores.

#2019 -------------------------------------------------------------------------------------------
As_DD_2019 <- read_delim("input/data/calidad del aire/2019/As_DD_2019.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

BaP_DD_2019 <- read_delim("input/data/calidad del aire/2019/BaP_DD_2019.csv", 
                          delim = ";", escape_double = FALSE, trim_ws = TRUE)

C6H6_DD_2019 <- read_delim("input/data/calidad del aire/2019/C6H6_DD_2019.csv", 
                           delim = ";", escape_double = FALSE, trim_ws = TRUE)

Cd_DD_2019 <- read_delim("input/data/calidad del aire/2019/Cd_DD_2019.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Ni_DD_2019 <- read_delim("input/data/calidad del aire/2019/Ni_DD_2019.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Pb_DD_2019 <- read_delim("input/data/calidad del aire/2019/Pb_DD_2019.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM10_DD_2019 <- read_delim("input/data/calidad del aire/2019/PM10_DD_2019.csv", 
                           delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM25_DD_2019 <- read_delim("input/data/calidad del aire/2019/PM25_DD_2019.csv", 
                           delim = ";", escape_double = FALSE, trim_ws = TRUE)

##Arsenico 2019

tablaAS_2019 <- select(.data = As_DD_2019, PROVINCIA, ANNO, D01:D31 )

tablaAS_2019 <-rename(.data = tablaAS_2019, AÑO = ANNO)

tabAS_2019 <- tablaAS_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabAS_2019$ARSENICO <- apply(tabAS_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabAS_2019 <- select(.data = tabAS_2019, PROVINCIA, AÑO, ARSENICO)

##Benzopireno 2019

tablaBaP_2019 <- select(.data = BaP_DD_2019, PROVINCIA, ANNO, D01:D31)

tablaBaP_2019 <-rename(.data = tablaBaP_2019, AÑO = ANNO)

tabBaP_2019 <- tablaBaP_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabBaP_2019$BENZOPIRENO <- apply(tabBaP_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabBaP_2019 <- select(.data = tabBaP_2019, PROVINCIA, AÑO, BENZOPIRENO)

##Benceno 2019

tablaC6H6_2019 <- select(.data = C6H6_DD_2019, PROVINCIA, ANNO, D01:D31)

tablaC6H6_2019 <-rename(.data = tablaC6H6_2019, AÑO = ANNO)

tabC6H6_2019 <- tablaC6H6_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabC6H6_2019$BENCENO <- apply(tabC6H6_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabC6H6_2019 <- select(.data = tabC6H6_2019, PROVINCIA, AÑO, BENCENO)

##Cadmio 2019

tablaCd_2019 <- select(.data = Cd_DD_2019, PROVINCIA, ANNO, D01:D31)

tablaCd_2019 <-rename(.data = tablaCd_2019, AÑO = ANNO)

tabCd_2019 <- tablaCd_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabCd_2019$CADMIO <- apply(tabCd_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabCd_2019 <- select(.data = tabCd_2019, PROVINCIA, AÑO, CADMIO)

##Niquel 2019

tablaNi_2019 <- select(.data = Ni_DD_2019, PROVINCIA, ANNO, D01:D31)

tablaNi_2019 <-rename(.data = tablaNi_2019, AÑO = ANNO)

tabNi_2019 <- tablaNi_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabNi_2019$NIQUEL <- apply(tabNi_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabNi_2019 <- select(.data = tabNi_2019, PROVINCIA, AÑO, NIQUEL)

##Plomo 2019

tablaPb_2019 <- select(.data = Pb_DD_2019, PROVINCIA, ANNO, D01:D31)

tablaPb_2019 <-rename(.data = tablaPb_2019, AÑO = ANNO)

tabPb_2019 <- tablaPb_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPb_2019$PLOMO <- apply(tabPb_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPb_2019 <- select(.data = tabPb_2019, PROVINCIA, AÑO, PLOMO)

##PM10 2019

tablaPM10_2019 <- select(.data = PM10_DD_2019, PROVINCIA, ANNO, D01:D31)

tablaPM10_2019 <-rename(.data = tablaPM10_2019, AÑO = ANNO)

tabPM10_2019 <- tablaPM10_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM10_2019$PM_10 <- apply(tabPM10_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM10_2019 <- select(.data = tabPM10_2019, PROVINCIA, AÑO, PM_10)

##PM25 2019

tablaPM25_2019 <- select(.data = PM25_DD_2019, PROVINCIA, ANNO, D01:D31)

tablaPM25_2019 <-rename(.data = tablaPM25_2019, AÑO = ANNO)

tabPM25_2019 <- tablaPM25_2019 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM25_2019$PM_25 <- apply(tabPM25_2019[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM25_2019 <- select(.data = tabPM25_2019, PROVINCIA, AÑO, PM_25)

#union de las 7 tablas 2019
join2019.1 <- full_join (x=tabAS_2019,y=tabBaP_2019, by=c("PROVINCIA","AÑO"))

join2019.2 <- full_join (x=join2019.1,y=tabCd_2019, by=c("PROVINCIA","AÑO"))

join2019.3 <- full_join (x=join2019.2,y=tabNi_2019, by=c("PROVINCIA","AÑO"))

join2019.4 <- full_join (x=join2019.3,y=tabPb_2019, by=c("PROVINCIA","AÑO"))

join2019.5 <- full_join (x=join2019.4,y=tabPM10_2019, by=c("PROVINCIA","AÑO"))

join2019.6 <- full_join (x=join2019.5,y=tabPM25_2019, by=c("PROVINCIA","AÑO"))


#2018 ------------------------------------------------------------------------------------- 

As_DD_2018 <- read_delim("input/data/calidad del aire/2018/As_DD_2018.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

B_a_P_DD_2018 <- read_delim("input/data/calidad del aire/2018/B(a)P_DD_2018.csv", 
                            delim = ";", escape_double = FALSE, trim_ws = TRUE)

Cd_DD_2018 <- read_delim("input/data/calidad del aire/2018/Cd_DD_2018.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Ni_DD_2018 <- read_delim("input/data/calidad del aire/2018/Ni_DD_2018.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Pb_DD_2018 <- read_delim("input/data/calidad del aire/2018/Pb_DD_2018.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM2_5_DD_2018 <- read_delim("input/data/calidad del aire/2018/PM2.5_DD_2018.csv", 
                            delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM10_DD_2018 <- read_delim("input/data/calidad del aire/2018/PM10_DD_2018.csv", 
                           delim = ";", escape_double = FALSE, trim_ws = TRUE)

##Arsenico 2018
tablaAS_2018 <- select(.data = As_DD_2018, PROVINCIA, ANNO, D01:D31 )

tablaAS_2018 <-rename(.data = tablaAS_2018, AÑO = ANNO)

tabAS_2018 <- tablaAS_2018 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabAS_2018$ARSENICO <- apply(tabAS_2018[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabAS_2018 <- select(.data = tabAS_2018, PROVINCIA, AÑO, ARSENICO)

##Benzopireno 2018
tablaBaP_2018 <- select(.data = B_a_P_DD_2018, PROVINCIA, ANNO, D01:D31)

tablaBaP_2018 <-rename(.data = tablaBaP_2018, AÑO = ANNO)

tabBaP_2018 <- tablaBaP_2018 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabBaP_2018$BENZOPIRENO <- apply(tabBaP_2018[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabBaP_2018 <- select(.data = tabBaP_2018, PROVINCIA, AÑO, BENZOPIRENO)

##Cadmio 2018

tablaCd_2018 <- select(.data = Cd_DD_2018, PROVINCIA, ANNO, D01:D31)

tablaCd_2018 <-rename(.data = tablaCd_2018, AÑO = ANNO)

tabCd_2018 <- tablaCd_2018 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabCd_2018$CADMIO <- apply(tabCd_2018[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabCd_2018 <- select(.data = tabCd_2018, PROVINCIA, AÑO, CADMIO)

##Niquel 2018

tablaNi_2018 <- select(.data = Ni_DD_2018, PROVINCIA, ANNO, D01:D31)

tablaNi_2018 <-rename(.data = tablaNi_2018, AÑO = ANNO)

tabNi_2018 <- tablaNi_2018 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabNi_2018$NIQUEL <- apply(tabNi_2018[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabNi_2018 <- select(.data = tabNi_2018, PROVINCIA, AÑO, NIQUEL)

##Plomo 2018

tablaPb_2018 <- select(.data = Pb_DD_2018, PROVINCIA, ANNO, D01:D31)

tablaPb_2018 <-rename(.data = tablaPb_2018, AÑO = ANNO)

tabPb_2018 <- tablaPb_2018 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPb_2018$PLOMO <- apply(tabPb_2018[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPb_2018 <- select(.data = tabPb_2018, PROVINCIA, AÑO, PLOMO)

##PM10 2018

tablaPM10_2018 <- select(.data = PM10_DD_2018, PROVINCIA, ANNO, D01:D31)

tablaPM10_2018 <-rename(.data = tablaPM10_2018, AÑO = ANNO)

tabPM10_2018 <- tablaPM10_2018 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM10_2018$PM_10 <- apply(tabPM10_2018[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM10_2018 <- select(.data = tabPM10_2018, PROVINCIA, AÑO, PM_10)

##PM25 2018

tablaPM25_2018 <- select(.data = PM2_5_DD_2018, PROVINCIA, ANNO, D01:D31)

tablaPM25_2018 <-rename(.data = tablaPM25_2018, AÑO = ANNO)

tabPM25_2018 <- tablaPM25_2018 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM25_2018$PM_25 <- apply(tabPM25_2018[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM25_2018 <- select(.data = tabPM25_2018, PROVINCIA, AÑO, PM_25)

#union de las 7 tablas 2018
join2018.1 <- full_join (x=tabAS_2018,y=tabBaP_2018, by=c("PROVINCIA","AÑO"))

join2018.2 <- full_join (x=join2018.1,y=tabCd_2018, by=c("PROVINCIA","AÑO"))

join2018.3 <- full_join (x=join2018.2,y=tabNi_2018, by=c("PROVINCIA","AÑO"))

join2018.4 <- full_join (x=join2018.3,y=tabPb_2018, by=c("PROVINCIA","AÑO"))

join2018.5 <- full_join (x=join2018.4,y=tabPM10_2018, by=c("PROVINCIA","AÑO"))

join2018.6 <- full_join (x=join2018.5,y=tabPM25_2018, by=c("PROVINCIA","AÑO"))


#2017 ----------------------------------------------------------------------------

As_DD_2017 <- read_delim("input/data/calidad del aire/2017/As_DD_2017.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

B_a_P_DD_2017 <- read_delim("input/data/calidad del aire/2017/B(a)P_DD_2017.csv", 
                            delim = ";", escape_double = FALSE, trim_ws = TRUE)

Cd_DD_2017 <- read_delim("input/data/calidad del aire/2017/Cd_DD_2017.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Ni_DD_2017 <- read_delim("input/data/calidad del aire/2017/Ni_DD_2017.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

Pb_DD_2017 <- read_delim("input/data/calidad del aire/2017/Pb_DD_2017.csv", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM2_5_DD_2017 <- read_delim("input/data/calidad del aire/2017/PM2.5_DD_2017.csv", 
                            delim = ";", escape_double = FALSE, trim_ws = TRUE)

PM10_DD_2017 <- read_delim("input/data/calidad del aire/2017/PM10_DD_2017.csv", 
                           delim = ";", escape_double = FALSE, trim_ws = TRUE)

##Arsenico 2017
tablaAS_2017 <- select(.data = As_DD_2017, PROVINCIA, ANNO, D01:D31 )

tablaAS_2017 <-rename(.data = tablaAS_2017, AÑO = ANNO)

tabAS_2017 <- tablaAS_2017 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabAS_2017$ARSENICO <- apply(tabAS_2017[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabAS_2017 <- select(.data = tabAS_2017, PROVINCIA, AÑO, ARSENICO)

##Benzopireno 2017

tablaBaP_2017 <- select(.data = B_a_P_DD_2017, PROVINCIA, ANNO, D01:D31)

tablaBaP_2017 <-rename(.data = tablaBaP_2017, AÑO = ANNO)

tabBaP_2017 <- tablaBaP_2017 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabBaP_2017$BENZOPIRENO <- apply(tabBaP_2017[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabBaP_2017 <- select(.data = tabBaP_2017, PROVINCIA, AÑO, BENZOPIRENO)

##Cadmio 2017

tablaCd_2017 <- select(.data = Cd_DD_2017, PROVINCIA, ANNO, D01:D31)

tablaCd_2017 <-rename(.data = tablaCd_2017, AÑO = ANNO)

tabCd_2017 <- tablaCd_2017 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabCd_2017$CADMIO <- apply(tabCd_2017[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabCd_2017 <- select(.data = tabCd_2017, PROVINCIA, AÑO, CADMIO)

##Niquel 2017

tablaNi_2017 <- select(.data = Ni_DD_2017, PROVINCIA, ANNO, D01:D31)

tablaNi_2017 <-rename(.data = tablaNi_2017, AÑO = ANNO)

tabNi_2017 <- tablaNi_2017 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabNi_2017$NIQUEL <- apply(tabNi_2017[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabNi_2017 <- select(.data = tabNi_2017, PROVINCIA, AÑO, NIQUEL)

##Plomo 2017

tablaPb_2017 <- select(.data = Pb_DD_2017, PROVINCIA, ANNO, D01:D31)

tablaPb_2017 <-rename(.data = tablaPb_2017, AÑO = ANNO)

tabPb_2017 <- tablaPb_2017 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPb_2017$PLOMO <- apply(tabPb_2017[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPb_2017 <- select(.data = tabPb_2017, PROVINCIA, AÑO, PLOMO)

##PM10 2017

tablaPM10_2017 <- select(.data = PM10_DD_2017, PROVINCIA, ANNO, D01:D31)

tablaPM10_2017 <-rename(.data = tablaPM10_2017, AÑO = ANNO)

tabPM10_2017 <- tablaPM10_2017 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM10_2017$PM_10 <- apply(tabPM10_2017[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM10_2017 <- select(.data = tabPM10_2017, PROVINCIA, AÑO, PM_10)

##PM25 2017

tablaPM25_2017 <- select(.data = PM2_5_DD_2017, PROVINCIA, ANNO, D01:D31)

tablaPM25_2017 <-rename(.data = tablaPM25_2017, AÑO = ANNO)

tabPM25_2017 <- tablaPM25_2017 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM25_2017$PM_25 <- apply(tabPM25_2017[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM25_2017 <- select(.data = tabPM25_2017, PROVINCIA, AÑO, PM_25)

#union de las 7 tablas 2017
join2017.1 <- full_join (x=tabAS_2017,y=tabBaP_2017, by=c("PROVINCIA","AÑO"))

join2017.2 <- full_join (x=join2017.1,y=tabCd_2017, by=c("PROVINCIA","AÑO"))

join2017.3 <- full_join (x=join2017.2,y=tabNi_2017, by=c("PROVINCIA","AÑO"))

join2017.4 <- full_join (x=join2017.3,y=tabPb_2017, by=c("PROVINCIA","AÑO"))

join2017.5 <- full_join (x=join2017.4,y=tabPM10_2017, by=c("PROVINCIA","AÑO"))

join2017.6 <- full_join (x=join2017.5,y=tabPM25_2017, by=c("PROVINCIA","AÑO"))


#2016
library(readxl)
As_DD_2016 <- read_excel("input/data/calidad del aire/2016/As_DD_2016.xlsx")

B_a_P_DD_2016 <- read_excel("input/data/calidad del aire/2016/B(a)P_DD_2016.xlsx")

Cd_DD_2016 <- read_excel("input/data/calidad del aire/2016/Cd_DD_2016.xlsx")

Ni_DD_2016 <- read_excel("input/data/calidad del aire/2016/Ni_DD_2016.xlsx")

Pb_DD_2016 <- read_excel("input/data/calidad del aire/2016/Pb_DD_2016.xlsx")

PM10_DD_2016 <- read_excel("input/data/calidad del aire/2016/PM10_DD_2016.xlsx")

PM25_DD_2016 <- read_excel("input/data/calidad del aire/2016/PM25_DD_2016.xlsx")

##Arsenico 2016

tablaAS_2016<- select(.data = As_DD_2016,PROVINCIA, ANNO, D01:D31 )

tablaAS_2016<-rename(.data = tablaAS_2016, AÑO = ANNO)

tabAS_2016<- tablaAS_2016 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabAS_2016$ARSENICO <- apply(tabAS_2016[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabAS_2016 <- select(.data = tabAS_2016, PROVINCIA, AÑO, ARSENICO)

##Benzopireno 2016

tablaBaP_2016<- select(.data = B_a_P_DD_2016,PROVINCIA, ANNO, D01:D31)

tablaBaP_2016 <-rename(.data = tablaBaP_2016, AÑO = ANNO)

tabBaP_2016<- tablaBaP_2016 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabBaP_2016$BENZOPIRENO <- apply(tabBaP_2016[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabBaP_2016 <- select(.data = tabBaP_2016, PROVINCIA, AÑO, BENZOPIRENO)

##Cadmio 2016

tablaCd_2016<- select(.data = Cd_DD_2016,PROVINCIA, ANNO, D01:D31)

tablaCd_2016 <-rename(.data = tablaCd_2016, AÑO = ANNO)

tabCd_2016<- tablaCd_2016 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabCd_2016$CADMIO <- apply(tabCd_2016[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabCd_2016 <- select(.data = tabCd_2016, PROVINCIA, AÑO, CADMIO)

##Niquel 2016

tablaNi_2016<- select(.data = Ni_DD_2016,PROVINCIA, ANNO, D01:D31)

tablaNi_2016 <-rename(.data = tablaNi_2016, AÑO = ANNO)

tabNi_2016<- tablaNi_2016 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabNi_2016$NIQUEL <- apply(tabNi_2016[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabNi_2016 <- select(.data = tabNi_2016, PROVINCIA, AÑO, NIQUEL)

##Plomo 2016

tablaPb_2016<- select(.data = Pb_DD_2016,PROVINCIA, ANNO, D01:D31)

tablaPb_2016 <-rename(.data = tablaPb_2016, AÑO = ANNO)

tabPb_2016<- tablaPb_2016 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPb_2016$PLOMO <- apply(tabPb_2016[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPb_2016 <- select(.data = tabPb_2016, PROVINCIA, AÑO, PLOMO)

##PM10 2016

tablaPM10_2016<- select(.data = PM10_DD_2016,PROVINCIA, ANNO, D01:D31)

tablaPM10_2016 <-rename(.data = tablaPM10_2016, AÑO = ANNO)

tabPM10_2016<- tablaPM10_2016 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM10_2016$PM_10 <- apply(tabPM10_2016[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM10_2016 <- select(.data = tabPM10_2016, PROVINCIA, AÑO, PM_10)

##PM25 2016

tablaPM25_2016<- select(.data = PM25_DD_2016,PROVINCIA, ANNO, D01:D31)

tablaPM25_2016 <-rename(.data = tablaPM25_2016, AÑO = ANNO)

tabPM25_2016<- tablaPM25_2016 %>%
  group_by(PROVINCIA,AÑO) %>%
  summarise(across(c(D01:D31), ~ mean(.x, na.rm = TRUE)))

tabPM25_2016$PM_25 <- apply(tabPM25_2016[ ,c(3:33)], 1, mean, na.rm = TRUE)

tabPM25_2016 <- select(.data = tabPM25_2016, PROVINCIA, AÑO, PM_25)

#union de las 7 tablas 2016
join2016.1 <- full_join (x=tabAS_2016,y=tabBaP_2016, by=c("PROVINCIA","AÑO"))

join2016.2 <- full_join (x=join2016.1,y=tabCd_2016, by=c("PROVINCIA","AÑO"))

join2016.3 <- full_join (x=join2016.2,y=tabNi_2016, by=c("PROVINCIA","AÑO"))

join2016.4 <- full_join (x=join2016.3,y=tabPb_2016, by=c("PROVINCIA","AÑO"))

join2016.5 <- full_join (x=join2016.4,y=tabPM10_2016, by=c("PROVINCIA","AÑO"))

join2016.6 <- full_join (x=join2016.5,y=tabPM25_2016, by=c("PROVINCIA","AÑO"))

Con todo el bloque de código anterior obtenemos 5 tablas. Cada una de ellas se corresponde a un año de calidad de aire con todas sus magnitudes. Ahora lo que tenemos que hacer es juntarlas.

Unión de todos los años en una única tabla.

Para ello, coo todas ellas tienen los mismos nombres en las columnas utilzaremos union_all.

#union de todas las tablas de años creadas anteriormente

Calidad1 <- union_all(join2020.6,join2019.6)

Calidad2 <- union_all(Calidad1,join2018.6)

Calidad3 <- union_all(Calidad2,join2017.6)

Calidad4 <- union_all(Calidad3,join2016.6)

DT::datatable(Calidad4)

A partir del código anterior obtenemos una única tabla con los datos de calidad de aire de los 5 años.

Cambio de la provincia a su nombre y comunidad.

En la columna provincia queremos cambiar el número que se corresponde a cada provincia, por el nombre de esta.

Para ello vamos a emplear una tabla con los dos datos, tanto el número como la provincia. Importamos la tabla con los números de provicia, su nombre y comunidad.

Hacemos un full_join a partir de la columna PROVINCIA.

library(readxl)
provincia_con_numero <- read_excel("input/data/calidad del aire/provincia con numero.xlsx", 
                                   sheet = "Estaciones evaluación 2019")

calidadFinal = full_join (x=Calidad4, y=provincia_con_numero, by=c("PROVINCIA"))

DT::datatable(calidadFinal)

Conseguimos una tabla calidadFinal. La cual es la tabla de Calidad4 anterior a la que se le han añadido 2 columnas correspondientes a los nombres de las provincias y a los de comunidades autónomas.

5.Datos de Enfermedades

Importación de los datos

Lo primero que realizamos es la importación de los datos de las enfermedades, lo hacemos a través de un archivo csv y un archivo txt. En realidad son los mismos sets de datos pero con formatos diferentes.

X49971.1 <- read_delim("input/data/Enfermedades/49971-2.csv", 
                       delim = ";", escape_double = FALSE, trim_ws = TRUE)
X49971.2 <- read_delim("input/data/Enfermedades/enf3.txt", 
                   delim = "\t", escape_double = FALSE, 
                   col_types = cols(Total = col_double(), 
                                    Total2 = col_double()), trim_ws = TRUE)

Unión de las partes que nos interesan.

El problema que nos ha surgido mientras realizabamos el trabajo es que la tabla X49971.1 considera los datos en formato character. Al pasar estos datos a numérico R considera el “.” como una coma y no como separador de miles, esto provoca que haya datos incorrectos en nuestro set.

Por ello, hemos decidido cambiar el formato de los datos directamente en un fichero excel sin puntos ni comas, se corresponde al fichero X49971.2. Con ello surge otro problema y es que los caracteres extraños como letras con tilde o ñ, aparecen con un símbolo como este �.

Hemos decidido por tanto, unir de cada una de las tablas una parte para conseguir una unica con los datos que nos interesan. Lo hacemos mediante el siguiente código:

TabEnf1 <- select(.data= X49971.1, Nacional:año)
TabEnf2 <- select(.data= X49971.2, Total)

#Union de ambas tablas
Enfermedades <- cbind(TabEnf1,TabEnf2)

Obtenemos una tabla con 7 columnas: Nacional, Comunidades y ciudades autónomas, Causa de Muerte, Sexo, Lugar de Ocurrencia, Año y Total. Obtemos un total de 214,200 entradas.

Esta tabla va a ser necesario modificarla para su posterior comparación. Lo que vamos hacer es cambiar el nombre de algunas columnas para evitar confusiones y que el acceso sea más facil.

Cambios en los nombres de columna.

colnames(Enfermedades) <- c('NACIONAL','CCAA','Causa_muerte','Sexo','Lugar','AÑO','Total')

Nombres correctos para comunidades.

Ahora vamos a cambiar el nombre de todas las comunidades a mayúscula para que coincida con los datos de calidad de aire, en realidad lo que hacemos es crear una nueva columna con los nombres en el formato correcto.

Enfermedades <- Enfermedades %>% 
  mutate(N_CCAA = case_when(CCAA == "Andalucía"  ~ "ANDALUCÍA",
                            CCAA == "Aragón"  ~ "ARAGÓN",
                            CCAA == "Asturias, Principado de" ~ "ASTURIAS (PRINCIPADO DE)",
                            CCAA == "Balears, Illes"  ~ "BALEARES (ISLAS)",
                            CCAA == "Canarias"  ~ "CANARIAS",
                            CCAA == "Cantabria"  ~ "CANTABRIA",
                            CCAA == "Castilla y León" ~ "CASTILLA Y LEÓN",
                            CCAA == "Castilla-La Mancha"  ~ "CASTILLA-LA MANCHA",
                            CCAA == "Cataluña"  ~ "CATALUÑA",
                            CCAA == "Comunitat Valenciana"  ~ "COMUNIDAD VALENCIANA",
                            CCAA == "Extremadura"  ~ "EXTREMADURA",
                            CCAA == "Galicia"  ~ "GALICIA",
                            CCAA == "Madrid, Comunidad de"  ~ "MADRID",
                            CCAA == "Murcia, Región de"  ~ "MURCIA (REGIÓN DE)",
                            CCAA == "Navarra, Comunidad Foral de"  ~ "NAVARRA (COMUNIDAD FORAL)",
                            CCAA == "País Vasco"  ~ "PAÍS VASCO",
                            CCAA == "Rioja, La"  ~ "RIOJA (LA)",
                            CCAA == "Ceuta"  ~ "CEUTA",
                            CCAA == "Melilla"  ~ "MELILLA"))

DT::datatable(Enfermedades)
## Warning in instance$preRenderHook(instance): It seems your data is too big
## for client-side DataTables. You may consider server-side processing: https://
## rstudio.github.io/DT/server.html

Eliminar NA en CCAA

Para que nuestras graficas en las preguntas posteriores aparezcan correctamente, tenemos que eliminar los NA en las comunidades para evitar que salgan los totales nacionales.

EnfermedadesFinal<-Enfermedades[!is.na(Enfermedades$N_CCAA),]

PREGUNTA 1: Efecto de la cantidad de PM 2,5 sobre las muertes por Neumonía.

Para realizar esta pregunta vamos a describir los conceptos:

  • Neumonía: Inflamación de los pulmones, causada por la infección de un virus o una bacteria, que se caracteriza por la presencia de fiebre alta, escalofríos, dolor intenso en el costado afectado del tórax, tos y expectoración.

  • PM2,5: Material particulado respirable presente en la atmósfera de nuestras ciudades en forma sólida o líquida (polvo, cenizas, hollín, partículas metálicas, cemento y polen, entre otras). Las PM2,5 de diámetro aerodinámico inferior o igual a los 2,5 micrómetros. Su origen está principalmente en fuentes de carácter antropogénico como las emisiones de los vehículos diesel.

Lo que vamos a buscar es si existe alguna relación entre las dos. Para ello, lo primero que vamos a hacer es seleccionar de nuestra tabla de EnfermedadesFinal la Neumonía.

EnfermedadesFinal %>%
  filter(Lugar == "Total") %>% 
  filter(Sexo == "Total") %>%
  filter(Causa_muerte=="063 Neumonía")-> Neumonia

DT::datatable(Neumonia)

Ahora vamos a hacer una media de nuestros datos de PM 2,5. Buscamos obtener un valor para cada comunidad y año.

estudioPM_25<- calidadFinal %>%
  group_by(N_CCAA,AÑO) %>%
  summarise(PM_25 = mean(PM_25, na.rm = TRUE))

DT::datatable(estudioPM_25)

Grafico de relación

ejercicio1 <- full_join (x=Neumonia, y= estudioPM_25,by=c("N_CCAA","AÑO"))

DT::datatable(ejercicio1)
##Gráfico
ggplot(data = ejercicio1, aes(x = PM_25, y = Total)) +
  geom_point(aes(colour = AÑO)) +
  stat_smooth() +
  theme_classic() +
  labs(
    x = "PM 2'5 (µg/m3)",
    y = "Número de muertes",
    title = 'Muertes por Neumonía vs PM 2,5',
    colour = 'Años'
  )

Explicación del gráfico 1:

Para responder a la primera pregunta sobre si la neumonía se ve afectada por los niveles de PM25 hemos reflejado esta grafica en la que compramos los niveles anuales de PM25 cada comunidad autónoma con el número de muertes anuales por neumonía de cada comunidad autónoma.

Una vez realizado el grafico pasamos a su análisis en el que podemos ver una clara correlación ya que cuanto mayor es la cantidad de PM25 hay más muertes por neumonía con esto podemos intuir que las pequeñas partículas sólidas o líquidas de polvo, cenizas, hollín, partículas metálicas, cemento o polen, dispersas en la atmósfera, y cuyo diámetro aerodinámico es menor que 2,5 µm tienen un efecto sobre la neumonía.

Un aspecto interesante a analizar es que parece como que al final se dispara la pendiente de la gráfica lo que podría indicar que hasta unos niveles como 12,5 los niveles de PM25 no son nocivos para el ser humano y no tiene influencia sobre la neumonía pero a partir de esos valores si que tienen un claro efecto negativo sobre la salud de las personas mas en concreto sobre la neumonía, aunque esto no se puede comprobar y para ello necesitaríamos de valores de una ciudad en la que los niveles de PM25 sean mas altos para intentar buscar esa correlación que hablamos.

PREGUNTA 2: Efecto de la cantidad de PM 10 sobre 4 enfermedades.

Para realizar esta pregunta vamos a describir los conceptos:

  • Las PM10 son partículas sólidas o líquidas de polvo, cenizas, hollín, partículas metálicas, cemento o polen, dispersas en la atmósfera, y cuyo diámetro varía entre 2,5 y 10 µm. Son compuestos inorgánicos como silicatos y aluminatos, metales pesados entre otros, y material orgánico asociado a partículas de carbono (hollín). Se caracterizan por poseer un pH básico debido a la combustión no controlada de materiales.

  • Asma: Es una enfermedad crónica que provoca que las vías respiratorias de los pulmones se hinchen y se estrechen. Esto hace que se presente dificultad para respirar como sibilancias, falta de aliento, opresión en el pecho y tos.

  • Neumonía: Inflamación de los pulmones, causada por la infección de un virus o una bacteria, que se caracteriza por la presencia de fiebre alta, escalofríos, dolor intenso en el costado afectado del tórax, tos y expectoración.

  • Bronquitis aguda y bronquiolitis: La bronquitis aguda es una infección a corto plazo de las vías respiratorias. La bronquiolitis es una infección de las vías respiratorias bajas que afecta a bebés y niños menores de 2 años. Es la causa más común de ingreso hospitalario de bebés menores de un año.

  • La gripe es una enfermedad infecciosa aguda, que afecta al aparato respiratorio y produce también una serie de síntomas generales característicos. Suele aparecer en brotes u oleadas, varias a lo largo de cada invierno, y es una enfermedad de distribución mundial.

Una vez definidos estos conceptos vamos a extraer las enfermedades anteriores de nuestra tabla de EnfermedadesFinal, filtrandolo para que coja el total de sexos y el total del lugar. Es decir, que los datos sean tanto hombres como mujeres y que la causa de muerte se haya producido en todos los lugares.

EnfermedadesFinal %>%
  filter(Lugar == "Total") %>% 
  filter(Sexo == "Total") %>%
  filter(Causa_muerte=="062 Influenza (gripe) (incluye gripe por virus de la influenza pandémica o zoonótica identificados)")-> Gripe

EnfermedadesFinal %>%
  filter(Lugar == "Total") %>% 
  filter(Sexo == "Total") %>%
  filter(Causa_muerte=="063 Neumonía")-> Neumonia

EnfermedadesFinal %>%
  filter(Lugar == "Total") %>% 
  filter(Sexo == "Total") %>%
  filter(Causa_muerte=="064 Enfermedades crónicas de las vías respiratorias inferiores (excepto asma)")-> EnfViasInf

EnfermedadesFinal %>%
  filter(Lugar == "Total") %>% 
  filter(Sexo == "Total") %>%
  filter(Causa_muerte=="065 Asma")-> Asma

Obtenemos así 4 tablas correspondientes a las 4 enfermedades a estudiar.

Ahora vamos a seleccionar de nuestra tabla de calidadFinal de aire las PM10, lo que necesitamos es tener un único valor por comunidad y por año. Creamos una tabla que se llame estudioPM_10 con la media de las PM1o agrupadas por año y comunidad.

estudioPM_10<- calidadFinal %>%
  group_by(N_CCAA,AÑO) %>%
  summarise(PM_10 = mean(PM_10, na.rm = TRUE))

Como queremos estudiar 4 enfermedades en un único grafico. Lo que vamos a hacer es un full_join de la tabla de cada enfermedad con la tabla de PM10. Una vez hecho esto unimos el resultado en una única tabla mediante union_all.

ejercicio2 <- full_join(x = estudioPM_10, y=Asma)
ejercicio2.1 <- full_join(x = estudioPM_10, y= Neumonia)
ejercicio2.2 <- full_join(x = estudioPM_10, y= Gripe)
ejercicio2.3 <- full_join(x = estudioPM_10, y= EnfViasInf)

ejercicio2a <- union_all(ejercicio2,ejercicio2.1)
ejercicio2b <- union_all(ejercicio2a,ejercicio2.2)
ejercicio2Fin <- union_all(ejercicio2b,ejercicio2.3)

DT::datatable(ejercicio2Fin)

Grafico de relación

ggplot(data = ejercicio2Fin, aes(x = PM_10, y = Total)) +
  geom_point() +
  stat_smooth() +
  theme_minimal() +
  facet_wrap( ~ Causa_muerte, nrow = 2)+
  labs(
    x = "PM 10 (µg/m3)",
    y = "Numero de muertes",
    title = 'Muertes por Enfermedades Respiratorias vs PM 10',
    colour = 'Años'
  )

Explicación del gráfico 2:

En este grafico creado con las partículas menores de 10µm contra el numero de muertes que ha habido en España durante 5 años por enfermedades que consideramos que afectan al sistema respiratorio.

Los resultados obtenidos podemos apreciar una correlación en la mayoría de ellas en la que a mayores cantidades de PM 10 aumentan las muertes por estas enfermedades salvo al final de estas cuando se alcanzan los niveles de PM10 mas altos que en este caso el numero de muertes se reduce en los 4 casos que hemos tomado.

Nuestra interpretación de estos gráficos es que hay una clara correlación entre las partículas del aire menores de 10µm y que el hecho de que esa correlación se pueda ver afectada cuando se alcanzan los valores máximos de PM10 y el numero de muertes baja se puede ver influenciado porque en los lugares donde hay peor calidad del aire (datos del aire con PM10 mas altos) hay menos habitantes por este factor de que la calidad del aire es baja y por consecuente hay un menor número de muertes.

Podemos concluir con que la exposición al PM10 puede provocar efectos nocivos en el sistema respiratorio

PREGUNTA 3: Efecto del Arsenico sobre la diabetes.

Queremos ver si el arsénico tiene algún efecto sobre la diabetes. Lo primero que vamos a hacer es definir los dos conceptos:

  • El arsénico es un elemento natural que se encuentra en la tierra y entre los minerales. Los componentes del arsénico se usan para preservar la madera, como plaguicidas y en ciertas industrias. El arsénico forma parte del aire, el agua y la tierra a través del polvo que se lleva el viento.

  • La diabetes es una enfermedad en la que los niveles de glucosa (azúcar) de la sangre están muy altos. La glucosa proviene de los alimentos que consume.

Ahora lo que vamos a hacer es poner un único valor de arsénico por cada comunidad y año:

estudioAs<- calidadFinal %>%
  group_by(N_CCAA,AÑO) %>%
  summarise(ARSENICO = mean(ARSENICO, na.rm = TRUE))

Tenemos que extraer de nuestra tabla de enfermedades los datos correspondientes a la diabetes, lo haremos mediante el siguiente bloque de código:

EnfermedadesFinal %>%
  filter(Lugar == "Total") %>%
  filter(Sexo == "Total") %>%
  filter(Causa_muerte=="044 Diabetes mellitus")-> Diabetes

Nos interesan los datos totales de lugar y total de sexo.

El siguiente paso es juntar los dos datos en una única tabla, para poder realizar el gráfico. Lo haremos mediante un full_join.

ejercicio3 <- full_join (x=Diabetes, y= estudioAs,by=c("N_CCAA","AÑO"))

Grafico de relación.

ggplot(data = ejercicio3, aes(x = ARSENICO, y = Total)) +
  geom_point() +
  stat_smooth() +
  theme_classic() +
  labs(
    x = "ARSENICO (µg/m3)",
    y = "Numero de muertes",
    title = 'Muertes por diabetes vs ARSENICO',
  )

PREGUNTA 4: Cancer vs Benzopireno y Arsénico

Lo que buscamos es ver cual de las dos magnitudes del aire tienen un efecto mayor sobre las muertes producidas por cáncer. Conceptos:

  • El arsénico es un elemento natural que se encuentra en la tierra y entre los minerales. Los componentes del arsénico se usan para preservar la madera, como plaguicidas y en ciertas industrias. El arsénico forma parte del aire, el agua y la tierra a través del polvo que se lleva el viento.

  • El benzopireno es un hidrocarburo policíclico aromático potencialmente carcinógeno (α-benzopireno) y que contienen algunos alimentos, como las carnes y el pescado.

  • El cáncer se caracterizan por el desarrollo de células anormales que se dividen sin control y tienen la capacidad de infiltrarse y destruir el tejido corporal normal.

Lo primero que vamos a hacer es coger los datos de Benzopireno y Arsénico de nuestra tabla de Calidad de Aire. Buscamos un único valor por comunidad y año.

estudioAs<- calidadFinal %>%
  group_by(N_CCAA,AÑO) %>%
  summarise(ARSENICO = mean(ARSENICO, na.rm = TRUE))

estudioBaP<- calidadFinal %>%
  group_by(N_CCAA,AÑO) %>%
  summarise(BENZOPIRENO = mean(BENZOPIRENO, na.rm = TRUE))

Ahora vamos a extraer de nuestra tabla de enfermedades los datos correspondientes a tumores.

EnfermedadesFinal %>%
  filter(Lugar == "Total") %>%
  filter(Sexo == "Total") %>%
  filter(Causa_muerte=="009-041 II.Tumores")-> Cancer

Ahora vamos a juntar las 3 variables a estudiar en una tabla para poder graficarlo

ejercicio4.a <- full_join (x=Cancer, y= estudioAs,by=c("N_CCAA","AÑO"))
ejercicio4 <- full_join (x=ejercicio4.a, y= estudioBaP,by=c("N_CCAA","AÑO"))

Grafico de relación.

ejercicio4 %>%
  select(Causa_muerte,Total,N_CCAA,ARSENICO,BENZOPIRENO)%>%
  pivot_longer(.,names_to = "Variable", values_to = "Valores", cols = c(ARSENICO:BENZOPIRENO)) -> ej4

DT::datatable(ej4)
ggplot(data = ej4, aes(x = Valores, y = Total)) +
  geom_point(aes(colour = factor(Variable)),na.rm = TRUE) +
  stat_smooth(na.rm = TRUE) +
  theme_light()+
  facet_wrap( ~ Variable, nrow = 2)+
  coord_cartesian(xlim = c(0, 0.73), expand = TRUE)+
  labs(
    x = "Valores (µg/m3)",
    y = "Numero de muertes",
    title = 'Muertes por cancer comparando niveles de ARSENICO Y BENZOPIRENO',
    colour = 'Variables'
  )

Explicación del gráfico 4:

En este estamos sacando la correlación que tiene el cáncer con algunas de los componentes del aire, en los que por su composición intuimos que pueden tener mayor relevancia en las muertes por cáncer en España.

Y claramente podemos observar como cada grafico tiene una tendencia diferente, en el caso de la relación de los niveles altos de arsénico con el número de muertes por cáncer, podemos apreciar una clara correlación en el que cuando se alcanzan los niveles más altos de arsénico también alcanzamos los niveles mas altos de muertes por cáncer en España, en cambio en el grafico construido con los datos de mortalidad por cáncer contra los niveles de benzopireno podemos observar que al inicio si que hay un crecimiento que puede llevar correlación entre el numero de muertes por cáncer y los niveles de benzopireno pero a medida que estos niveles de benzopireno crecen más el numero de muertes desciendes levemente.